rm(list=ls())

(WD <- getwd())
if (!is.null(WD)) setwd(WD)

indata0      = paste0(WD,"/_other_data/_nz/","", collapse = NULL)

library(readxl)
library(foreign)
library(haven)
library(AER)
library(sandwich)
library(lmtest)
library(pracma)
library(dplyr)
library(stargazer)
library(plm)
library(pracma)
library(DataCombine)
library(jtools)
library(plyr)
library(ggplot2)

lagpad <- function(x, k) {
  res <- c(rep(NA, k), x)[1:length(x)]
  return(res)
}

unc         = 1.11 # forecaster uncertainty in std. 
unc1        = 0.89 # forecaster uncertainty in std others.
dis         = 3.06 # disagreement in std.
dis1        = 2.43 # disagreement in std others.

pub_unc     = 1.68 # public forecaster uncertainty in std

data        = read_excel(paste(indata0,"new_zealand.xls",sep=""), sheet = "Data")
tmp          = row(data)
tmp          = tmp[,1]
data$qdate   = 1961.00+(tmp-1)/4
data         = data[data$qdate>=1970.00,]
data         = data[data$qdate<=2020.00,]

data_short   = data[data$qdate>=1985.00,]
var_inf      = var(data_short$inflation)
var_inf_use  = var_inf

w_hat_private_1      = unc1/unc
w_hat_private_2      = dis1/dis

var_inf                  = var_inf_use

precision_private_others = w_hat_private_1/(1-w_hat_private_1)*(1/var_inf)
uncertainty_others_1      = sqrt(1/(precision_private_others+var_inf^(-1)))

precision_private_others = w_hat_private_2/(1-w_hat_private_2)*(1/var_inf)
uncertainty_others_2      = sqrt(1/(precision_private_others+var_inf^(-1)))

uncertainty = rep(NA, 2)
uncertainty[1] = unc
uncertainty[2] = unc1

stdev       = rep(NA, 2)
stdev[1]    = dis
stdev[2]    = dis1

implied_unc_1 = rep(NA, 2)
implied_unc_1[1] = uncertainty_others_1

implied_unc_2 = rep(NA, 2)
implied_unc_2[1] = uncertainty_others_2

table_3_all  =  data.frame(uncertainty, stdev, implied_unc_1,implied_unc_2)
col         = c('Inflation one-year ahead','Consensus inflation expectations')
table_3_all$col = col
write.table(format(table_3_all,digits=3), file = "_tables/table_2_bottom.txt", sep = "\t")


